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ABSTRACT 


The  nonlinear  stress-strain  behavior  of  unidirectional  fiber 
composite  material  was  examined.  Two  approaches  were  adopted.  The 
first  approach  considers  the  stress-strain  as  nonlinear  elastic  relations, 
which  are  derived  at  through  the  use  of  a  complementary  energy  density 
function,  and  takes  into  account  the  material  symmetry.  The  second 
approach  adopts  the  Ramberg-Osgood  representation  of  one-dimensional 
stress-strain  curves.  By  utilizing  the  deformation  theory  to  express  the 
nonlinear  parts  of  the  strain,  a  nonlinear  stress-strain  behavior  was 
obtained. 

An  incremental  formulation  of  these  theories  was  implemented  in  a 
finite  element  analysis.  These  approaches  were  examined  in  the  context  of 
postbuckling  analysis  of  thick-section  laminated  plates.  It  is  shown  that 
the  nonlinear  behavior  leads  to  a  stiffness  reduction  for  thick-section 
plates.  Moreover,  the  models  examined  produce  a  close  response  for  thick 
multi-layered  laminated  plates. 
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1.  INTRODUCTION 


The  design  of  fiber  composite  laminated  structures  has  grown 
rapidly  over  the  last  two  decades.  These  laminated  anisotropic 
structures  offer  an  exquisite  alternative  to  the  conventional  isotropic 
materials  due  to  their  light  weight  to  high  strength  ratio.  Therefore, 
it  is  not  surprising  to  find  that  the  most  immediate  and  primary 
applications  of  laminated  structures  were  exclusive  to  the 
aeronautical  and  aviation  industry  in  the  form  of  thin  laminated 
plates  and  shells.  However,  recently  there  has  been  a  growing 
interest  in  the  design  and  application  of  fiber-reinforced  composite 
materials  in  the  form  of  thick  section  plate  and  shell  structures.  This 
introduces  a  major  concern  when  using  methods  of  analysis,  such  as 
the  classical  laminate  theory  or  its  modifications,  since  these  methods 
assume  the  presence  of  a  plane  stress  field  in  each  lamina.  However, 
in  a  thick  section  composite,  we  verge  to  a  true  tridimensional  stress 
state,  especially  if  a  need  to  analyze  stress  concentrations  or 
structural  discontinuities  arises. 

The  search  for  effective  and  accurate  design  methods  to 
analyze  laminated  structures  has  led  to  the  use  of  complex  models. 
These  complexities  arise  through  the  consideration  of  geometric  and 
material  nonlinear  behavior.  Such  models  are  assumed  to  produce  a 
reliable  and  accurate  design. 
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Recently,  more  attention  has  been  directed  towards  exploring 
the  effects  of  nonlinear  constitutive  behavior  in  the  application  of 
thermoset,  thermoplastic  and  metal  matrix  fiber  composite  systems. 
These  systems  demand  the  application  of  a  proper  nonlinear 
constitutive  theory.  A  nonlinear  constitutive  model  can  have  a 
significant  effect  on  the  structural  design,  especially  in  the  presence 
of  local  stress  concentrations,  such  as  crack  tips,  holes  and  cutouts, 
edge  effects  and  structural  discontinuities.  Hence,  it  is  important  to 
assess  such  effects  with  the  proper  choice  of  constitutive  relations. 

This  study  is  concerned  with  the  formulation  of  a  nonlinear 

constitutive  model  for  fiber  composite  laminates.  The  nonlinear 
behavior  of  a  single  lamina,  which  results  from  the  matrix  nonlinear 

behavior  or  from  a  presence  of  micro-effects  such  as  micro-cracks 

and  imperfect  fiber  matrix  interface,  is  introduced  in  the  lamina's 

constitutive  model. 

In  the  first  part  of  this  report,  two  major  theories  and 
approaches  of  nonlinear  constitutive  models  are  examined  and 
reviewed  in  a  detailed  manner.  In  the  second  part,  an  incremental 
formulation  of  the  above  theories  is  implemented.  A  nonlinear  finite 
element  analysis  is  performed  to  examine  the  above  theories  in  the 
context  of  postbuckling  analysis  of  thick  section  laminated  plates. 
Discussion  and  recommendations  follow. 
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2.  LITERATURE  REVIEW 


2.1  Nonlinear  Behavior  of  Fiber  Composite  Laminates 

It  has  been  well  recognized  that  different  fiber  composite 
systems  exhibit  some  form  of  nonlinear  stress-strain  behavior.  One 
source  of  nonlinearity  can  be  related  to  the  matrix  material,  such  as 
viscoelastic  behavior  in  polymeric  matrices  or  plasticity  in  metal 
matrices.  Another  source  of  nonlinearity  can  result  due  to  micro- 
structural  characteristics,  such  as  micro-cracks  or  imperfect  interface 
between  the  fiber  and  the  matrix. 

Several  levels  and  approaches  can  be  used  in  modeling  the 
nonlinear  behavior  of  fiber  composites.  A  micromechanics  approach 
can  be  obtained  by  modeling  the  micro-structure  of  the  fibers  and 
matrix  inclusion  in  some  form  of  periodical  array  subject  to  different 
states  of  stress  and  strain,  Adams  [1].  Another  approach  is  to 
consider  the  nonlinear  responses  of  a  homogenized  anisotropic 
lamina  and  to  examine  the  nonlinear  responses  for  the  one 
dimensional  average  transverse  stress  and/or  for  the  axial  shear, 
while  the  stress-strain  relation  in  the  fiber  direction  remains  linear. 
Finally,  interaction  terms  of  multiaxial  stress  states  can  be  formed 
through  the  use  of  a  potential  stress  function  to  describe  the 
nonlinear  parts  of  the  strains,  such  as  the  stress  energy  density 
function  or  the  deformation  theory  of  plasticity. 
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IIachin  et  al  [2]  formulated  and  used  this  approach  whereby 
they  characterized  the  nonlinear  stress-strain  relations  in  a 
Ramberg-Osgood  representation  for  the  one  dimensional  average 
transverse  and  axial  shear  stress-strain  relations.  Moreover,  they 
used  the  deformation  theory  of  plasticity  to  form  nonlinear 
constitutive  relations  which  included  an  interaction  term  between 
nonlinear  axial  shear  and  transverse  stresses. 

Smith  and  Rivlin  [3j  considered  a  nonlinear  polynomial  form  of 
the  stress  or  strain  energy  density  function,  in  which  they  reduced 
the  polynomial  basis  by  employing  symmetry  conditions  on  the  basis 
form  and  by  using  some  mathematical  group  theorems.  It  should  be 
noted  that  this  approach  included  major  crystal  classes  of  anisotropy. 

The  key  idea  in  applying  the  above  approach  to  anisotropic 
materials,  is  to  extend  the  symmetry  that  a  material  exhibits  in  its 
crystallographic  form  to  anisotropic  bodies  that  are  not  crystals,  but 
have  a  micro-structural  symmetry,  which  results  from  the  fiber 
periodic  positioning  in  its  inclusion,  disregards  the  matrix  or  fiber 
materials. 

Hahn  and  Tsai  [4]  followed  this  approach  and  considered  the 
axial  shear  nonlinear  behavior  in  a  composite  lamina  under  a  plane 
stress  state.  They  used  a  complementary  elastic  energy  density 
function  and  imposed  transverse-isotropy  symmetry  conditions  to 
generate  a  nonlinear  elastic  constitutive  relation.  The  interaction 
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terms  between  the  shear  nonlinear  stress  and  other  stress 

components  were  neglected  in  this  formulation. 

The  theories  of  Hashin  et  al  and  Smith  and  Rivlin  are  reviewed 
and  implemented  in  an  incremental  nonlinear  finite  element  analysis 
in  the  following  chapters. 

2.2  Postbuckling  Analysis  of  Laminated  Plates 

A  fiat  plate  subject  to  membrane  forces  remains  flat  provided 
these  forces  are  small.  By  increasing  the  applied  forces  the  flat  plate 
becomes  unstable  and  tends  to  deform  in  such  configuration  to  allow 
small  lateral  displacements.  These  loads  are  the  critical  or  buckling 
loads.  If  the  loads  are  increased  more,  the  plate  establishes  a  stable 
configuration  which  allows  bending  due  to  membrane  forces.  The 

plate  is  then  said  to  be  in  the  postbuckling  range  and  a  stiffness 
reduction  will  occur,  but  the  plate  can  still  resist  the  increasing 
applied  loads. 

2.2.1  Nonlinear  Analysis  of  Thin  Laminated  Plates 

The  involvement  of  geometric  nonlinearity  along  with 
laminated  plate  theory  makes  it  a  tedious,  if  not  impossible,  task  to 
generate  close  form  solutions  in  the  postbuckling  range.  The 

postbuckling  formulations  of  thin  flat  plates  involve  the  use  of 

equilibrium  equations  written  in  the  current  deformed  configuration, 
along  with  some  Airy  type  function  generating  equilibrium 
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conditions  in  terms  of  the  membrane  forces.  A  second  nonlinear 
equation  is  derived  from  the  compatibility  condition  to  form  a  set  of 
two  coupled  nonlinear  partial  differential  equations  (FDE).  This 
formulation  of  governing  postbuckling  equations  for  thin  flat  plates 
is  introduced  in  Appendix  B.  These  nonlinear  PDE's  can  be  solved 
approximately  through  the  use  of  a  truncated  double  Courier  series  to 
generate  a  set  of  nonlinear  algebric  equations. 

This  approach  was  used  by  [5]  to  solve  a  postbuckling  of 
orthotropic  rectangular  plate  with  all-clamped  or  all-simply 
supported  edges.  Also  [  6  j  considered  a  postbuckling  of 
symmetrically  laminated  plates  with  the  same  boundary  conditions 
as  1 5]. 


2.2.2  Nonlinear  Analysis  of  Thick  Composite  Laminate 

The  analysis  of  thick  section  laminated  plates  and  shells 
includes  the  account  for  transverse  shear  which  becomes  an 
important  factor  as  the  depth  to  length  ratio  increases.  The 
transverse  shear  included  by  the  Reissner-Mindlin  plate  theory 
offers  a  better  alternative  to  the  Kirchhoff  plate  theory.  Postbuckling 
analysis  of  thick  composite  laminates  was  considered  by  17,8].  Both 
Kirchhoff  thin  plate  and  Mindlin  thick  plate  were  analyzed.  It 

indicates  that  the  transverse  shear  reduces  the  buckling  loads, 
moreover,  the  postbuckling  stiffness  is  also  reduced  by  including  the 

transverse  shear  effects.  The  major  assumption,  in  the  nonlinear 

range,  of  both  Kirchhoff  and  Mindlin  type  plate  theories  is  the 
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The 


existence  of  small  strains  and  large  in-plane  rotations, 
independent  interpolation  of  the  rotations,  in  the  x  and  y  directions, 
irom  the  transverse  deflection  in  Mindlin  plate  theory  makes  it  a 

difficult  task  to  formulate  the  governing  postbuckling  equations.  The 
finite  element  method  of  Mindlin  type  elements  offer  an  attractive 

method  for  the  solution  of  postbucklmg  analysis  of  laminated  plates. 

Finally,  it  is  important  to  notice  that  various  combined 

experimental  and  numerical  studies  were  performed  to  validate  the 
thin  laminated  plate  theory  in  uniaxial  compression  over  the 

postbuckling  range  [9],  To  the  knowledge  of  the  authors,  no  biaxial 
compressive  experimental  and  numerical  study  was  conducted  to 
verify  the  theory. 
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3.  NONLINEAR  CONSTITUTIVE  MODELS  OF  FIBER  COMPOSITE 

LAMINATES 


The  unusual  shear  behavior  in  thermoset  fiber  composites  and 

matrix  nonlinearity  in  thermoplastic  composite  systems  led  to  a 
concern  to  assess  the  nonlinear  behavior  through  a  suitable 

formulation  of  constitutive  model  and  through  the  examination  of 
such  nonlinear  effects  on  the  structural  responses. 

In  this  chapter,  two  major  theories  are  reviewed  and  examined. 
The  first  theory  considers  the  stress-strain  behavior  as  a  nonlinear 
elastic.  A  complementary  energy  density  function  is  formed  and 
various  classes  of  anisotropic  symmetry  are  imposed  on  its 
polynomial  basis.  The  result  is  a  simplified  reduced  polynomial 
basis,  which  a  nonlinear  constitutive  model  can  be  derived  from.  The 
second  approach  used  the  deformation  theory  of  plasticity  and  the 

Ramberg-Osgood  characterization  of  a  nonlinear  stress-strain 
relation.  These  constitutive  models  are  both  implemented  for  the 

AS4/J  polymer  composite  material. 

3.1  Anisotropic  Elastic  Constitutive  Model 

3.1.1  Strain  Energy  and  Complementary  Energy  Density  Functions 

When  given  a  perfectly  elastic  body,  there  exists  an  elastic 
strain  energy  function  per  unit  volume,  W,  of  the  undeformed  body. 
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This  function  is  a  single-valued  function  of  the  strain  tensor  such 
that 


1  id w  aw\  aw 

Oj ;  =  —  I - +  —  Oj,  =  — 

2  \3eij  dcji)  or  aeij  (1) 

where  W  =  W(e;j)  =  W(uij)  (2) 

The  complementary  energy  density  function  W*  can  be  defined 
through 


W(ejj)  +  W*(Oij)  =  Ojj  Cjj  (  3 ) 

It  is  important  to  notice  that  under  symmetry  operations  Ojj  ejj 
remains  invariant.  Hence,  one  can  choose  a  representation  of  the 
constitutive  model  through  W  or  W*. 

3.1.2  Polynomial  Representation  of  Strain  Energy  Function 

Assuming  that  the  strain  energy  function  is  analytic  in  its 
arguments,  it  can  be  represented  as 


W  (ejj)  =  X  Wr(eij) 

r=i  (4) 


Where  Wr  is  a  homogeneous  function  of  the  degree  r 
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W^e.j)  =  XrWr(eij) 


(5) 


Obviously  it  is  impossible  to  consider  infinite  terms  of  Wr. 
Hence,  a  truncated  cubic  approximation  will  be  used.  It  can  be 
written  as 


W  —  W„  +  Cjj  ejj  +  ^  Cjjk]  Cjj  ek]+  ^  jK)mn  eij  eki  emn 

Where  W0  is  an  arbitrary  constant  and  Cjj  is  an  initial  stress 
tensor.  It  is  easy  to  notice  that  a  quadratic  form  of  W  yields  linear 
constitutive  relations. 

3.1.3  Material  Symmetry 

If  the  anisotropic  material  is  considered  to  be  homogeneous 
and  possesses  symmetry  due  to  the  fibers  positioning,  certain 
restrictions  can  be  imposed  on  the  polynomial  basis  of  W.  Given  that 
X;  and  Xj  are  two  equivalent  symmetric  orientations,  one  can  write 


W(eij)  =  WCe'ij)  =  W(aik  aji  ekl) 


(7) 


and  aij  =  cos(x'i,  Xj) 

Next,  we  will  define  some  basic  symmetric  transformation 


1  0 


operators 


I  --  Identity  transformation. 


-I  --  Inverse  transformation. 

Rj  ..  Reflection  of  the  i-th  axis,  xi  =  -  xi. 

Sln,a  Right-handed  rule  rotation  through  a  degrees  around  the 
axis  in  the  direction  of  the  normal  n. 

These  symmetric  operators  will  be  used  to  define  some  basic 
crystalline  systems. 

Triclinic  System  (I,  -I) 

The  identity  and  inverse  transformations  do  not  effect  the 
stress  or  strain  tensors.  Hence,  there  can  be  no  restrictions  placed  on 
the  polynomial  basis  of  the  strain  energy  function. 


W(ejj)  =  W(en,  e22,  633,  e^,  e23>  en) 


(8) 


Monoclinic  Symmetry 

For  material  having  Monoclinic  symmetry,  the  transformation 
operators  defining  this  symmetry  are 

I,  -I,  S'3,71  (9) 

The  first  two  operators  do  not  effect  the  form  of  the  polynomial 
basis.  The  rotation  of  the  cartesian  coordinate  system  in  71  degrees 
around  the  X3  axis  yields 


^23  =  'e23 


and 


®13  =  'e23 


(10) 


Hence,  the  .  nitation  imposed  on  the  form  of  the  strain  energy 
function  is 

W(en,  e22,  e33,  ei2,  e23,  ei3>  =  W(en,  e22,  e33,  -e23.  -en)  (11) 

This  will  create  two  sets  of  variables  which  the  polynomial 
basis  is  symmetric  in 

(ai,  a2)  =  (e23,  ei3)  and  (bi,  b2)  =  (-e23,  -en)  (12) 

Applying  group  theorems  (1)  and  (3)  in  appendix  A  to  produce 
a  polynomial  basis  for  W  in  the  form 

W  =  W Ou  h,  h,  U,  Ki,  K2)  Kn,  K12,  K22)  (13) 

Where  Jj  are  the  invariant  quantities  under  the  above 
transformation,  and 

Kj  =  K2  =  0 

Kn=e23,  Ki2=ei3C23,  ^23=^23 

The  reduced  polynomial  basis  of  the  strain  energy  is 
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(14) 


W  =  W(en,  Z22 ,  e33,  c\2,  e]3,  e^3>  eo  e23) 


(15) 


Orthotropic  Symmetry 

A  material  having  three  mutually  orthogonal  planes  of 
symmetry  is  called  orthotropic  material.  The  basic  symmetry 
operators  are 


I,  -I,  s^'n,  Si2-n,  si3'n 


(16) 


The  first  three  operators  were  applied  to  produce  Equation 
(15),  applying  the  fourth  operator  yields 


e12  =  -en 


e23  =  -e23 


(17) 


The  symmetry  limitation  imposed  on  W  in  equation  (15)  will 


be 


W(en,  e22.  e33,  e^,  e23,  ef3,  ei3  e23)  = 


W(en,  e22,  e33,  -ej2,  e|3>  ef3,  -en  e23) 


(18) 


We  can  identify  two  sets  of  variables  which  the  polynomial 
basis  is  symmetric  in 

Cai ,  a2)  =  (ei2,  en  e23)  and  (bi>  62)  =  (-ei2,  -en  e23) 
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(19) 


Using  group  theorems  (1)  and  (3)  in  appendix  A,  the 
polynomial  basis  of  W  is 


W  =  W(Jl5  J2,  J3,  U,  Is.  K!,  K2,  K„,  K12,  K22) 


(20) 


Where  J[  are  the  transformation  invariants,  writing  Kj  and  Kjj 


Kj  =  K2  =  0 


Ki2  =  ei2ei3e23? 


K22  =  ej3e23  (21) 


Hence, 


W  =  W(en,  e22,  e33,  ej2,  e23,  e^3,  ei2ei3e23) 


(22) 


Applying  S'1-71  will  not  effect  the  form  of  W  because  if  there  are 
two  orthogonal  symmetric  planes  the  third  orthogonal  plane  is  also 
symmetric. 

Transverse  Isotropy 

A  composite  unidirectional  fibrous  reinforced  material  aligned 
in  x3  direction,  in  which  the  Xj  x2  plane  is  compsoed  of  a  hexagonal  or 
random  array  of  fibers,  is  said  to  have  a  transverse  isotropic 
symmetry.  The  basic  symmetric  operators  are 


I,  -I,  Ri,  Sl3« 


(a  -  arbitrary) 


(23) 


For  an  arbitrary  rotation  of  the  coordinate  system  around  x3 
axis,  the  new  coordinate  system  can  be  expressed  as 


X;  =  a;j  Xj 


(24) 


An  alternative  representation  of  (24)  is 


X{  +iX^  =  e-ia(X!  +iX2), 


xj=x3 


(25) 


Can  show 


ell  +  e22 


en  +  e22 


e33 


e33 


ell  ■  e22  +  2ie12 


e-2ia  (en  -  e22  +  2i  ei2) 


(26) 


e13  +  ,e23 


e-ia  (ei3  +  ie23) 


We  can  identify  the  following  invariants 


Ji=en+e22  J2-e3i+e23 


J3  =  e33 


(27) 


If  we  define  the  quantities 
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Ei  =  en  -  e22 
E2  =  2en 


Fi  =  ej3  -  e^3 
F2  =  2  en  e23 


(28) 


From  the  transformation,  can  show 

Et  =  Ej  cos  2a  +  E2  sin  2a 

E2  =  E2  cos  2a  -  Ei  sin  2a 

(29) 

Fj  =Fi  cos  2a  +  F2  sin  2a 

F2  =  F2  cos  2a  -  Fi  sin  2a 

Symmetry  requirement  imposed  on  the  form  of  W,  can  be 
stated  as 

W(Ji,  J2,  J3,  Ei,  E2,  Fi,  F2)  =  W(Ji,  J2,  J3,  El5  E2,  Fj,  F2)  (30) 

Applying  group  theorem  (2)  in  Appendix  A,  where  the  vectors 

a(1)  ja(2) 

wi  and^i  are 

a;1 2  =  (Ei,  E2)  ap}  -  (Fj,  F2)  n  =  2  (31) 

The  invariant  forms  due  to  the  above  orthogonal 
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transformation  are 


+  Ej  =  Ji  -4J4 


aP)ap)  =  E1Fi+E2F2  =  Jl 


(32) 


ap)ap)  =  Fl2  +  Fj  =  2J5  +  ixJ2  -  2J3J4 


A  =  EiF2  -  E2F] 


Where  J 1  ^  J2,  J3  aredefined  in  (27)  and 


J4  =  en  e22  - 1\2 


J5  —  lejjl 


(33) 


W  can  be  expressed  in  the  polynomial  basis  as 


W  =W(Ji,  J2,  J3,  J4,  J5) 


(34) 


Assumptions  and  Simplifications 

Consider  a  transversely  isotropic  material  where  Xi  coincides 
with  the  fibers  direction.  If  the  given  lamina  is  to  be  considered  in  a 
state  of  plane  stress,  the  form  of  W  will  be  considerably  simplified  as 


W  =  W(en,  e22,  e]2) 
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(35) 


The  same  argument  can  be  made  with  the  complementary 
stress  energy  density  function,  using  the  relation  (3) 

W*  =  W*(Oi,  02,  al)  (36) 

An  approximated  fourth  degree  of  W*  can  be  written  as 
W*  =  i  Snof  + 1  S22  °2  +  S12  ai  o2  +  S66  ^6  + 

jSmo]  +^-S222  ^2  +  Si12  a?  o2  +  Sm  Oi  <*2  + 

S  1 66  1  0^  +  S266  02<56  +  ^Siin  oj  +  ^  S2222^2  + 

"  S6666  <*6  +  Sui2  Oi  <*2  +  S1122  <5\  ©I  + 

Sl222  Ol  <*2  +  S1 166  O?  ©I  +  S2266  <*2  <*6  +  S1266  <*1  C2  eg 

According  to  Hahn  and  Tsai  [4],  the  interaction  terms  between 
a,  and  o6  or  o2  and  o6  are  neglected  for  the  material  they 
experimented  with.  It  has  been  shown  that  an  important  effect  is 
produced  on  the  resolved  transverse  stress  in  a  uniaxial  test  for  a 
wide  range  of  fiber  orientations.  Hence,  the  only  effect  considered 
by  their  model  was  the  effect  of  shear  stress  nonlinearity  without 
any  interaction  terms  which  may  involve  other  stress  components  in 
the  form  of  W*. 
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Finally,  the  complementary  energy  density  function  is  reduced 


to 


W*=£S„  c?i  +  S22  02  +  S12  Oi  o2  +  2  S66  o^  +  Is6666ot 


'38) 


Note  that  the  shear  components  appeal  separately  in  their 
powers  and  that  there  are  no  interaction  terms.  Also,  the  only 
nonlinearity  is  in  the  axial  shear.  The  constitutive  stress-strain 
relations  can  be  derived  from 


eij  “ 


aw* 
a  o;, 


(39) 


3.2  Deformation  Theory  and  Ramberg-Osgood  Representation 

Deformation  type  nonlinear  stress-strain  relations  can  be 
represented  in  the  form 


®ij  —  Sjjki  CTy  (40) 

Where  Sijk|  is  a  function  of  the  average  stresses  or  the  proper 
stress  invariants.  If  the  material  is  transversely  isotropic,  and  the 
fibers  are  alligned  in  the  x,  direction,  the  stress  invariants  are 
obtained  in  Equations  (27)  and  (33).  It  is  convenient  to  express  the 
strain  tensor  into  a  linear  elastic  part  e?j  and  a  nonlinear  strain  part  en 


(41) 


eij=eij+eij 


it  follows 


Sjjkl 


(42) 


where 


Sfjkl  =  SPjkl(J|,J2,J3,J4,J5) 


(43) 


Consider  a  thin  lamina  in  a  state  of  plane  stress,  the  effective 
nonlinear  compliances  are  a  function  of  the  three  average  plane- 
stress  components 


^ijkl  Sijk](t5 11  Gil®  i  2^ 

The  nonlinear  stress-strain  relations  are 

en  =S|iOn  +  S,2  022 

e22  =  °H  +  S22°22 


(44) 


(45) 


ei2  =  O12 


For  unidirectional  fiberous  reinforced  materials,  the  fibers 
direction  is  clearly  stiffer  than  the  matrix.  Therefore,  the  average 
stress  °n  can  be  assumed  to  be  carried  by  the  fibers.  The  other 
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stress  comnonents  G 22  and  Gi2  are  carried  by  the  matrix  and  give  the 
nonlinear  behavior.  The  only  nonlinear  parts  of  the  compliance 
tensor  are 


^22  “  S22(022>G?2)  anc|  Sg6  -  S£6(022>of2) 


(46) 


Next,  one  dimensional  representation  of  nonlinear  stress-strain 
relations  in  Ramberg-Osgood  form  are 


C22 


e12 


1  + 


O^tN-l 
Ty  / 


(47) 


Where  E2  is  the  transverse  elastic  Young’s  modulus  and  G]2is 
the  axial  elastic  shear  modulus.  The  parameters  °y  Ty,  M,N  are  curve 
fitting  parameters. 


3.2.1  J2  -  Deformation  Theory 


The  basic  assumption  of  the  deformation  theory  of  plasticity, 
for  isotropic  materials,  is  that  the  plastic  strains  can  be  represented 

..  s 


efj  =  f(J2)  Sij 


(48) 


Where  S,,  is  the  deviatoric  stress  and  J2  is  its  second  invariant 


Sij  —  Oij  -  ^ 


(49) 


h  =  ^  S,jS'j 

Similarly,  we  can  generalize  this  to  anisotropic  materials  by 
assuming  that 


eij  =  f  (L)  Sij 


(50) 


Where  L  is  a  general  quadratic  function  of  the  stress  invariants. 
Consider  a  state  of  plane  stress,  the  invariant  forms  resulting  from 
symmetry  consideration  are 


Oii,o22,o?2  (51) 

Assuming  that  the  stress  an  does  not  produce  any  plastic 
strain  since  it  is  mainly  carried  by  the  fibers,  L  is  a  quadratic 
function  of  o22  and  °i2 


L  =  L  (a22.o]2)  =  a2<*22  +  P2°?2  (52) 

Note  that  there  is  nothing  fundamental  in  the  assumption  of 
the  general  quadratic  form  in  Equation  (50).  Next,  we  can  determine 
the  coefficients  a,  (3  from  a  comparison  of  the  nonlinear  one 
dimensional  stress-strain  relations  with  the  Ramberg-Osgood  stress- 
strain  representation. 
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A  Ramberg-Osgood  representation  of  one  dimensional 


transverse  stress 


e22  “  $22(022)  022  -  ^ 


022 

Gy 


M-l 


(53) 


or 


~p  022 

22  “eT  l 


(54) 


From  the  deformation  theory,  the  nonlinear  part  of  the 
transverse  strain  is 


*=^22(0=!^  felony 

fc-2  tL2 


Comparing  Equations  (54)  and  (55)  leads  to 


a2  =  -i-  and  f22  =  L*!1 


o  y 


(55) 


(56) 


The  function  f22  is  the  quadratic  form  L  to  the  (M-l)/2  power. 
Similarly,  in  the  case  of  one  dimensional  o12,  the  nonlinear 
strain  in  the  Ramberg-Osgood  form  is 


(57) 


From  the  deformation  theory,  the  nonlinear  shear  strain  is 
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(58) 


Pp  _ 
e12  “ 


=  ^f66(I.)  =  g!if66(|32 


'12 


'12 


Comparing  Equations  (57)  and  (58),  it  follows  that 

P2  =  \  and  f66  =  Lt 
*y 


(59) 


Finally,  the  nonlinear  compliances  can  be  expressed  as 


M-l 


S22  =  J“[(«^22)2+  (P<Ti2)2  ]  2 

c-2 


N-l 


S66  =  pM  {aa22?  +  (pcfn)2  }  2 
U12 


(60) 


w  here 


a 


2  —  _L 


2__L 

t2 
ty 


The  parameters  oY,  M  and  ty,N  respectively  are  curve  fitting 
parameters  for  the  Ramberg-Osgood  representation  of  stress-strain 
curves  in  a  one  dimensional  case. 

It  follows  that  the  nonlinear  plane  stress-strain  relations  are 


4.  FINITE  ELEMENT  FORMULATIONS 


The  objective  of  this  chapter  is  to  describe  the  finite  element 
method  used  to  solve  nonlinear  structural  systems.  The  well  known 
virtual  work  principal  is  used  as  a  'weak  form*  to  derive  the 
equilibrium  equations 


J  Bt  o  dv  -  P  =  0 


(62) 


Where  P  is  the  applied  force  vector,  B  is  the  strain- 
displacement  matrix  assembled  in  the  global  structural  coordinates 
along  with  the  stress  vector.  Equation  (62)  can  be  written  as 

I(u)  =  P  (63) 

Where  I  is  the  internal  resisting  force  vector  and  is  a  nonlinear 
function  of  the  nodal  displacements  in  the  case  of  nonlinear 
geometric  problems.  The  representation  of  the  virtual  work 

expression  in  Equation  (63)  can  lead  to  two  basic  finite  element 
formulations.  The  first  is  the  secant  stiffness  matrix  formulation 
which  relates  the  total  displacement  with  the  total  load  vector.  The 
second  method  is  an  incremental  formulation  which  defines  a 
relation  between  increments  of  loads  and  displacements  through  a 
tangent  stiffness  matrix.  This  method  can  be  derived  by 
differentiating  Equation  (63),  it  follows  that 


Kt  AU  =  Ap 


(64) 


where 


Kt  =  Ko+(K!  +Kj+  K2)  +  KG 
K0  =  |v  BfCBjdv 

=  L  B7CBnidv  (65) 

K2  =  jv  B^CBn,dv 
Kg  =  Jy  GTMcGdv 

Note  that  the  matrix  K,  is  nonsymmetric  and  is  a  linear  function 
in  the  nodal  displacements.  K0is  a  constant  independent  matrix.  The 
quantity  in  parenthesis  is  called  the  initial  displacement  matrix  and 
is  symmetric.  Different  combined  incremental-iterative  techniques 
are  used  to  solve  the  nonlinear  set  of  algebraic  equations  (64).  The 
Newton-Raphson  (NR)  and  its  modifications  are  the  most  known  and 
used  methods.  These  methods  can  be  described  as 

Kt  AU"1!  =  Pn+i  - 1„+1  (66) 

Where  i  is  the  iteration  number  and  n  is  the  step  or  increment 
number.  The  internal  resistance  vector  In+i  is  updated  according  to 
the  method  used. 
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The  load  increment  Pn+i  is  specified  for  a  fixed  load  level 
method,  whereas,  in  the  case  of  variable  NR  iterative  methods,  the 
load  increment  Pn+]  is  part  of  the  solution.  These  methods  follow  the 
structure  equilibrium  path  of  a  structure  initially  subject  to  a  given 
load  pattern  P0  in  the  form 


P  =  ^P0  (67) 

were  X  is  a  scalar  determined  by  adding  a  constraint  to  the 
equilibrium  equations. 

In  order  to  accurately  simulate  the  decrease  in  the  load  and 
displacement  near  limit  or  bifurcation  points,  where  the  tangent 
stiffness  becomes  singular,  a  modified-Riks  (11)  iterative  algorithm 
is  used  to  follow  the  equilibrium  path.  The  basic  idea  of  the  Riks 
method  is  moving  along  the  tangent  line  to  the  equilibrium  path  and 
searching  for  an  equilibrium  point  located  in  the  plane  orthogonal  to 
the  tangent  line. 

This  study  makes  use  of  the  ABAQUS  finite  element  code  (12) 
where  the  modified  Riks  method  is  implemented.  The  nonlinear 
formation  is  the  updated  Lagrangian  method  in  which  the  reference 
configuration  is  updated  to  an  intermediate  configuration  between 
the  original  and  current  configurations.  Thus,  the  current  and 
reference  configurations  remain  close. 
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4.1.  Element  Formulation 


The  plate  element  used  in  the  postbuckling  analysis  is  the 
isoparametric  8-nodes  quadrilateral  plate  element.  It  has  six 
degrees  of  freedom  (D.O.F.)  in  each  node  which  includes  three 
displacements  and  three  rotations.  The  element  is  based  on  the 
Reissner-Mindlin  plate  theory.  Numerical  formulations  are  more 
simple  and  convenient  to  implement  than  the  Kirchhoff  plate 
elements.  The  theory  also  accommodates  for  transverse  shear 
strains.  The  continuity  requirement  is  C°  continuity  due  to  the  fact 
that  only  first  order  derivatives  of  the  D.O.F.  appear  in  the  virtual 
work  expressions.  The  displacements  and  rotations  are  decoupled. 
Hence,  the  element  is  mappable  and  allows  modeling  of  complex 
geometries. 

The  displacement  field  is 

u(x,y,z)  =  u°(x,y)  -  Z0X 

v(x,y,z)  =  v°(x,y)  -  Z0y  (68) 


w(x,y,z)  =  w(x,y} 

The  above  displacement  field  is  independently  interpolated 
within  the  element 


u°  =  NjUi  V°  =  NjVj  W  =  NjWj 

(69) 

9X  =  Nj0*j  0y  =  NjOyj 
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The  notation  used  is  the  index  summation  on  the  element 


nodes.  The  nonlinear  strain-displacement  relations,  defined  in 
appendix  B,  are  applied  to  the  displacement  field  in  (68)  and  yield 


£xx  -  u°x  +  2  w?x  ■  Z  0X 


eyy  =  v?y  +  2  w?y  '  Z  0y 

£Xy  —  U‘.y  +  V^x  +  W,xW,y  -  Z  0x,y  -  Z0y  x 
^xz  =  w,x  -  0x  y 

£yZ  —  W,y  -  0  y 


(70) 


Both  Equations  (69)  and  (70)  defines  the  strain-displacement 
matrix  B  in  the  form 


8e*x 

"  Ni,x 

0 

0 
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0 

Sfiy  y 
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Equation  (71)  defines  linear  and  nonlinear  parts  of  the  strain- 
displacement  matrix.  Again,  note  the  summation  over  indices  i  and  j 
where  they  represent  the  element  nodes. 

4.2  Constitutive  Tangent  Stiffness 


The  material  tangent  stiffness  or  Jacobian  matrix  of  the 
constitutive  model  relates  incremental  stresses  and  strains.  It  can  be 
defined  as 


SAG; 

CiJ  =  6Aq 

Where  Acjj  and  Aej  are  the  stress  and  strain  incremental  vectors. 
In  the  case  of  material  and  geometric  nonlinearities,  the  material 
model  should  be  formed  in  an  incremental  form  to  accommodate 
with  the  incremental  finite  element  formulations. 


(72) 


Three  constitutive  models  are  considered  in  this  study.  The 
first  model,  presented  in  Chapter  3,  considers  the  axial  shear 
nonlinearity.  The  incremental  form  of  this  model,  in  a  plane-stress 
condition,  is 


S 12  0 

do. 

S22  0 

< 

do2  > 

0  S  66  +  3  S  666  6^6  - 

-do6. 

(73) 
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The  second  model  considers  uncoupled  transverse  and  shear  stress- 
strain  relations.  Its  incremental  form  is 


'den 
'  de2  r 
-de6J 


Si  i 

S  1  2 

1 

0 

'do,l 

II 

S  12 

S22  +  3  S  2222^2  0 

< 

do2 1 

-  0 

0 

S66  +  3  S  6666^6  - 

.do6J 

(74) 


The  third  model  considered  in  3.1.4  for  coupled  transverse  and  shear 
stress-strain  in  Ramberg-Osgood  representation.  The  incremental 
form  for  this  model,  in  the  case  of  M=N=3,  is 


'de,' 

<  de2 

11 

kde6- 

1 

S„  S12  0 

S , 2  S2i1+3a2a2  +  P2^)  2S22p2a2a6 

0  2S66a2o2o6  S66(l  +  a2°2  +  P2  ol)  J 


do, 

da2 

daj 


(75) 


The  Jacobian  matrix  of  this  model  is  not  symmetric.  This  will 
invoke  unsymmetric  equation  solves.  All  these  models  are  in  the 
materials  coordinate  system  and  need  to  be  transformed  in  the 
global  structural  coordinates  for  each  layer.  The  relation  for  the 
transverse  stress-strain  is  considered  linear  and  introduced  to  the 
formulation  through  shear  correction  factors. 

The  material  considered  in  this  study  is  AS4/Ji  polymer 
composite.  The  stress-strain  curves  are  presented  in  Fig.  2  and  Fig. 
3. 
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4.3  Imperfection  Analysis 


In  order  to  generate  a  solution  for  plates  in  the  postbuckling 
range  and  to  enable  the  plate  to  deform  in  the  transverse  direction 
due  to  inplane  forces,  an  imperfection  technique  is  implemented  in 
the  finite  element  model.  This  method  is  carried  out  by  either  a 
perturbation  of  the  model  nodes  in  the  transverse  direction  such  that 
it  will  include  some  of  the  first  eigen  vectors  in  its  form,  or  by 
applying  small  distributed  loads  in  the  transverse  direction  at  the 
beginning  of  the  analysis. 
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5.  NUMERICAL  RESULTS 


5.1  Plate  Geometry  and  Boundary  Conditions 

This  study  considers  the  problem  of  buckling  and  postbuckling 
analysis  of  laminated  plates.  The  material  constitutive  models 
developed  in  Chapter  3  are  implemented  in  an  incremental  nonlinear 
finite  element  formulation.  The  plates  have  a  square  geometry  (L  = 
300mm)  and  are  subjected  to  a  compressive  biaxial  in-plane  loading. 
The  laminate  stack  sequence  considered  is  a  quasi-isotropic  laminate 
[0/90/±45 |ns.  The  plate  edges  are  constrained  to  produce  a  uniform 
displacement  in  the  x  and  y  direction  respectively.  Both  in  the 
buckling  and  postbuckling  analysis  a  full  finite  element  mesh  is  used 
due  to  the  nonlinear  geometry  involved.  Most  of  the  results  arc 
presented  in  the  form  of  the  applied  load,  normalized  by  the  plate 
critical  load,  versus  end-shortening  in  the  x  direction,  normalized  by 
the  plate  thickness.  The  effect  of  stress-biaxiality  ratio  is  studied  for 
a  clamped  plate  with  different  aspect-ratios  (t/L  =  0.02,  0.04,  0.08). 
The  range  of  biaxial  loading  ratios  are  chosen  for  the  values:  0,  0.5, 
l.,2.  .  For  Ny/Nx=  0.  this  is  the  case  of  uniaxial  loading. 

5.2  Linearized  Buckling  and  Imperfection  Analysis 

In  order  to  enable  the  plate  to  deform  in  the  lateral  direction 

due  to  membrane  forces,  a  geometric  imperfection  is  introduced  in 

the  finite  element  mesh.  The  first  three  eigenmodes,  derived  from 
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linearized  buckling  analysis,  are  combined  in  a  linear  combination 
and  scaled  by  some  percentage  of  the  plate  thickness.  It  is  a  well 
known  fact  that  plate  structures  are  not  sensitive  to  the  form  of 
imperfection.  Nevertheless,  a  smaller  percentage  of  imperfection 
implemented,  the  more  the  postbuckling  response  is  close  to  the  real 
behavior.  One  percent  of  the  thickness  is  chosen  to  scale  the  linear 
combination  of  the  first  three  eigenmodes.  This  will  enable  the  plate 
response  to  be  very  close  to  the  real  perfect  structure. 


5.3  Materials  Models 

The  fiber  reinforced  material  considered  in  this  study  is  made 
of  AS4/J1  thermoplastic  matrix.  The  effective  elastic  properties  of 
this  composite  system  are: 

En  =  17.9  x  106  psi,  E22  =  0.9  x  106  psi,  i312  =  0.313 
Gl2  =  GI3  =  0.77  x  106  psi,  G23  =  0.31  x  106  psi 

The  nonlinear  stress-strain  curves  for  this  material,  in  axial 
shear  [13)  and  transverse  tensile-compression,  are  presented  in  figs. 
2-3.  The  three  nonlinear  constitutive  models  developed  in  Chapter  3 
are  formulated  based  on  these  one-dimensional  stress-strain  curves. 
These  models  are  referred  to  as  model  1,  2  and  3  in  the  following 
discursions. 
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5.4  Convergence  Study 


The  convergence  rate,  in  finite  element  analysis  of  plates, 
depends  strongly  upon  the  plate  aspect-ratio  (t/L).  A  large  number 
of  elements  are  required  to  obtain  an  accurate  solution  for  thick 
laminates,  especially  in  the  postbuckling  range.  Hence,  the  case 
study  chosen  is  a  clamped  plate  under  biaxial  compression  (Ny/Nx  = 
0.5)  with  an  aspect-ratio  of  0.08.  Fig.  4  illustrates  the  rate  of 
convergence  in  the  prebuckling  and  postbuckling  ranges  for  n  x  n 
element  mesh  refinement.  Based  on  these  results,  the  8x8  elements 
mesh  was  chosen  in  the  following  analysis.  It  should  be  noted  at  this 
stage  that  there  is  an  enormous  amount  of  computational  time 
needed  for  the  postbuckling  analysis  where  the  material  model  is 
nonlinear.  A  typical  8x8  isoparametric  8-nodes  quadrilateral  shell 
elements  modeling  a  [0/90/±45  ] 5s  laminated  plate  require  a  large 
number  of  Gaussian  integration  points  (2  x  2  x  64  x  40  Gaussian 
points).  This  calculation  assumes  that  each  lamina  is  modeled 
through  the  thickness  with  one  integration  point.  The  material 
model  would  be  called  about  10,000  times  per  one  iteration!  Ni  Ong 
that  a  typical  postbuckling  analysis  requires  hundreds  of  iterations, 
one  can  clearly  realize  the  large  amount  of  computational  time 
involved.  Hence,  the  stack  sequence  chosen  for  most  of  this  study 
was  |0/90/±45 1,.  A  design  of  thick  laminates  requires  the  use  of  a 
sublaminate  technique  whereby  a  typical  repeating  sublaminate 
includes  the  material  nonlinearities. 


5.5  Effect  of  Thickness 


The  laminate  aspect-ratio  effect  on  the  postbuckling  behavior 
is  studied  for  (0/90/±45]ns  clamped  plate  under  biaxial  compression 
(Ny/Nx  =  0,  0.5,.l.,  2.).  The  thickness/length  ratios  considered  are: 
0.02,  0.04,  0.08.  Figs.  5-8  illustrates  the  end-shortening  curves  of 
three  laminates  with  different  aspect-ratio  and  different  biaxial 
loading.  It  is  shown  that  the  postbuckling  stiffness  decreases,  as  the 
load  increases,  for  large  aspect-ratio  (t/L  =  0.08).  This  is  due  to  the 
growing  effect  of  transverse  shear  in  thick  plates.  The  first  nonlinear 
material  model  is  plotted  along  each  curve.  The  nonlinear  curves  of 
plate  response  with  aspect-ratio  less  than  0.04  do  not  show  a 
difference  compared  to  the  linear  material  response.  Note  the  early 
separation  of  the  linear  and  nonlinear  curves  for  t/L  =  0.08.  This 

indicates  that  there  is  a  difference  in  the  buckling  load  between  the 

two  models.  This  shift  or  separation  of  curves  shows  that  the 
nonlinear  plate  is  more  flexible.  This  flexibility  is  expressed  by 
looking  at  a  fixed  load  level  and  considering  that  the  corresponding 

difference  of  the  displacement  is  much  higher  than  the  difference  in 

the  load  for  a  fixed  value  of  the  displacement. 

It  can  be  concluded  that  the  more  the  postbuckling  stiffness  is 
decreasing,  the  more  the  effect  of  material  nonlinearities  can  be  seen. 
Also,  the  nonlinear  response  is  associated  with  large  deflections.  The 
only  nonlinear  material  model  considered  in  the  above  figures  is  the 
first  model  which  includes  the  axial  shear  nonlinearity.  The  other 
models  are  considered  in  the  following  figures  for  t/L  =  0.08. 
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5.6  Thick-Section  Plates 


The  three  nonlinear  models  are  considered  in  the  analysis  of  a 
thick-section  plate  (t/L  =  0.08)  with  biaxial  loading  (Ny/Nx  =  0/  0.5,  1, 
2)  and  stack  sequence  [0/90/±45]s.  The  end-shortening  of  the  plate 
versus  the  normalized  applied  load  is  presented  in  figs.  9-11.  The 
plate  center  deflection  is  presented  in  figs.  13  ’5.  These  figures 
clearly  demonstrate  the  effect  of  nonlinear  behavior  of  fiber 
composite  materials  upon  the  global  response  of  the  structure.  The 
closest  nonlinear  behavior  to  the  linear  response  is  the  nonlinear 
axial  shear  model  (model  1).  Model  2  and  3  follow  in  order.  This 
consistent  pattern  of  behavior  can  be  explained  by  the  fact  that 
model  2  includes  the  transverse  nonlinear  response  along  with  the 
nonlinear  in-plane  shear,  in  an  independent  fashion.  Model  3 
considers  both  nonlinearities  along  with  interaction  terms.  It  can  be 
expected  that  model  3  produces  the  largest  nonlinear  effect.  The 
same  conclusion  can  be  drawn  from  the  lateral  deflection  curves  figs. 
13-15  for  different  biaxial  loading  ratios. 

Finally,  an  additional  case  of  multi-layered  thick  laminate 
[0/90/±45  ]  12s  is  considered.  The  aspect  ratio  in  (t/L  =  0.08)  and  the 
biaxial  loading  ratio  is  Ny/Nx  =  1.  .  The  end-shortening  postbuckling 

response  is  presented  in  fig.  16.  the  lateral  displacement  behavior  is 
presented  in  fig.  17.  Both  curves  indicate  that  the  three  nonlinear 
models  produce  a  very  close  response  in  most  of  the  postbuckling 

range.  This  behavior  can  be  explained  by  the  presence  of  large 
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number  of  layers  of  0°  and  90°  degree  fiber  angles.  These  layers 
carrv  most  of  the  applied  loads  in  the  fiber  directions.  As  a  result, 
the  nonlinear  strain  distribution  in  the  matrix  is  reduced.  Figs.  18- 
19  support  this  observation.  In  these  figures  the  transverse  shear- 
stress  resultants  Vx  and  Vy  are  presented  respectively.  The  location 
of  these  resultants  is  the  mid-point  of  the  upper  right  quarter  of  the 
plate.  These  figures  indicate  that  the  nonlinear  models  do  not  affect 
the  stress  distribution  in  the  structure  as  much  as  producing 
relatively  large  displacements.  Also,  fig.  20  shows  the  strain  in  the 
fiber  direction  of  the  upper  and  lower  layers  at  the  plate  center.  It  is 
shown  that  the  strain  in  the  fiber  direction  can  be  considered  as  a 
linear  elastic  strain  with  no  plastic  or  nonlinear  part. 
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6.  CONCLUDING  REMARKS 


The  approach  utilized  in  this  report  was  to  present  two  major 
theories  of  nonlinear  constitutive  models  of  stress-strain  behavior. 
The  first  approach  considers  the  stress-strain  as  nonlinear  elastic 
relations,  which  are  derived  through  the  use  of  a  complementary 
energy  density  function  and  takes  into  account  the  material 
symmetries.  The  second  approach  was  to  adopt  the  Ramberg-Osgood 
representation  of  the  one-dimensional  nonlinear  stress-strain  curves. 
By  utilizing  the  deformation  theory  as  a  representation  of  the 
nonlinear  parts  of  the  strain,  a  nonlinear  stress-strain  behavior  was 
obtained.  The  first  two  models,  referred  as  model  1  and  model  2,  are 
based  on  the  first  approach.  Model  1  represents  the  nonlinear  axial 
shear  behavior.  In  addition,  model  2  takes  into  account  the 
transverse  stress-strain  behavior  in  an  independent  fashion.  The 
third  model  is  based  on  the  second  approach  of  deformation  theory. 
It  includes  both  nonlinear  stress-strain  relations.  Moreover,  an 
interaction  expression  was  formulated  to  account  for  simultaneous 
presence  of  axial  shear  and  transverse  stress. 

The  postbuckling  results  of  thick  section  laminated  plates 
indicates  that  the  nonlinear  stress-strain  behavior  has  a  global  effect 
on  the  structural  response.  The  stiffness  was  reduced  and  a  large 
displacement  was  associated  with  the  structural  behavior.  The 
distribution  of  the  stress  resultant  in  plate  seems  not  to  have  a  major 
change  compared  to  the  linear  response. 


It  is  important  to  emphasize  that  the  problems  analyzed  in  this 
study  do  not  include  any  stiess  singularities  such  as  matrix  cracking, 
delamination,  holes  and  cutouts.  These  local  effects  can  have  a 
significant  inelastic  behavior.  The  question  of  the  effect  of  nonlinear 
constitutive  behavior  on  such  local  damage  modes  and  structural 
discontinuities  should  be  addressed  as  a  continuation  of  this  study. 

A  comparison  of  the  three  models  presented  yields  that  the 
first  model  produces  the  most  conservative  response.  The 
interaction  terms  in  the  deformation  theory  model  gives  a  difference 
of  response  than  the  other  models.  The  structural  behavior  of  multi¬ 
layered  thick  laminates  which  are  associated  with  these  constitutive 
models  are  very  close.  This  result  indicates  that  the  presence  of 
large  numbers  of  unidirectional  layers  oriented  in  different 
directions  decreases  the  nonlinear  strains  by  the  fact  that  the  fiber 
direction  can  carry  most  of  the  applied  loads. 

Failure  criterias  were  not  used  in  this  study.  Obviously,  the 
laminated  thick  section  plate  will  be  damaged  in  the  process  of 
postbuckling  or  even  prebuckling  response.  Hence,  it  is  important  to 
include  failure  and  local  effects  in  a  future  study. 
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APPENDIX  A 


The  following  short  review  covers  three  mathematical  group 
theorems.  Mathematical  discussion  and  proofs  are  beyond  the 
context  of  this  study  and  can  be  found  in  [10]. 

Theorem  (1) 

A  polynomial  basis  for  polynomials,  which  are  symmetric  in 
the  two  sets  of  variables  (a,,  a2,  ....,  an)  and  (bi,  b2,  ....  bn),  is  formed 
from  the  quantities 

Kj  =  ^(aj+bj) 

j.k  =  1,2 . n 

Kjk  =  2-(ajbk  +  akbj) 


Theorem  (2) 

A  polynomial  basis  for  polynomials,  which  are  symmetric  in 
the  n  vectors  a(r)  =  («i,(X2,---,ccn)(r)  in  n-dimensional  space,  and  which  is 
invariant  under  all  proper  orthogonal  transformations,  is  formed  by 

the  scalar  products 

a(r)a(s) 

and  i,r,s  =  1,2,  ...,n 

ki 

Theorem  (3) 


A  polynomial  basis  for  a  polynomial  in  the  variables  (ai,a2,...,an), 
(b!,b2,...,bn)  and  (Ji,J2,...,Jk),  which  is  form-invariant  under  a  group  of 
transformations  and  (J,,J2,...,Jk)  are  invariants,  is  formed  by  adjoining 
to  the  quantities  (Ji,J2,...,Jk)  the  polynomials  in  the  variables 
(a!,a2,...,an)  and  (bi,b2,...,bn)  which  are  form-invariant  under  the  given 
group  of  transformations. 
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APPENDIX  B 


Postbucking  Analysis  of  Thin  Plates 

The  purpose  of  the  present  appendix  is  to  develop  the 
governing  equations  for  postbuckling  analysis  of  laminated  thin 
plates.  Although  in  this  study  the  postbuckling  analysis  was 
performed  through  the  use  of  the  finite  element  method  and  for 
thick  and  thin  plates,  it  is  important  to  recognize  the  formulation  of 
the  governing  equations  for  postbuckling  analysis  and  that  there  are 
different  ways  to  obtain  approximated  solutions,  such  as  the  finite 
difference  method  or  by  a  truncated  double  fourier  series  method. 

In  attempting  to  write  the  equilibrium  equations  for 
postbuckling  of  plates,  one  should  write  the  equilibrium  equations  in 
the  current  configuration.  By  transforming  these  equations  to  the 
original  or  reference  configuration,  a  set  of  nonlinear  equilibrium 
equations  can  be  obtained. 

A  point  occupying  the  reference  configuration  is  denoted  as  a 
position  vector  x,  this  vector  corresponds  to  the  position  vector  x  in 
the  deformed  configuration,  where 

Xi  =  Xj  +  u;  ( 1 ) 


Where  Uj  is  the  displacement  vector, 
tensor  is  defined  as 

r-  3Xj  r. 

F>j=^-  =  5ij  +  Ui.j 


The  deformation  gradient 


(2) 
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The  Green  strain  tensor  can  be  written  as 


eij  =  ^{FkiFkj  -  5ij)  =  +  uj.i  +  Uk.jujc.j)  3  j 

Equilibrium  in  the  deformed  configuration,  with  no  body  forces,  can 
be  expressed  as 

(4) 

Where  Ti  and  ds  are  traction  vector  and  element  area  respectively  in 
the  deformed  configuration.  The  stress  measure  used  is  the 
symmetric  second  Piola-Kirchhoff  stress  which  is  obtained  by 
mapping  the  force  transmitted  over  an  element  area  ds  in  the  current 
configuration  to  the  reference  configuration 

Ti  ds  =  Fy  Tj  ds  =  Fjj  ajk  nk  ds  (5 ) 


Where  Tj  is  a  fictitious  stress  vector  acting  over  an  element  area  ds, 
with  a  unit  outward  normal  nk,  in  the  original  configuration. 
Substituting  Equation  (5)  into  (4)  yields 


{ Fjj  ajk  ds  =  0 


Using  the  divergence  theorem 

}y(FijOjk),k  dv  =  0 


(6) 


(7) 
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Equation  (7)  is  true  for  any  volume  v.  The  equilibrium  equations  are 


substituting 


[FijOjklk  —  [(§ij  +ui,j  )  t?jk],k  —  0 


(8) 


i  =  1,2,3 


Uij  =  ejj  +  Wij 


(9) 


Where  e,j  is  the  linear  engineering  strain  and  is  the  skew- 
symmetric  rotation-tensor 


[(8ij  +  eij  +  w;j)ojk],k  -  0  (10) 

Next,  we  will  use  the  kinimatic  assumptions  of  nonlinear  deformation 
of  thin  plates,  which  are  known  as  the  von  Karman  theory  of  plates. 
This  theory  assumes  that  the  strains  can  be  neglected  compared  with 
the  rotations.  Moreover,  for  flat  thin  plates  we  can  neglect  the 
rotation  w12.  Introducing  these  assumptions  in  Equation  (10) 

[(8ij  +  wjj  )cjjk],k  =  0  i  =  1,2,3  (11) 

Writing  the  set  of  Equations  (11)  explicitly,  we  can  show  that  the 
first  tow  equations,  i  =  1,2,  are  reduced  to  the  well  known  linear 
equilibrium  equations.  The  only  nonlinear  equation  is  the 
equilibrium  equation  in  the  Z  direction 

aap,p  =  0  a,P  =  1,2  (12) 


The  third  equilibrium  equation  is 


(13) 


(w,x^xx  +  w,y^yx  +  ^zx),x  -f  (wiX(5Xy  +  w,y^yy  +  ^zy),y 

4-  (w,x^xz  "*■  w,yOyz  +  O'zzJ.z  =  0 


Where  w  is  the  lateral  displacement.  Integrating  the  above 
equations  in  the  same  fashion  as  when  developing  the  classical 


Kirchhoff  plate  equations 


N„aa  —0 


(14) 


Mxx.xx  +  2MXyiXy  +  Myyyy  +  NXX  Wxx  +  2NXy  Wjy 
+  Nyy  Wyy  +  P  =  0 


(15) 


The  nonlinear  strains  are 


£xx  =  ufx  +  w?x  +  Zkxx 

£yy  =  V^y  +  ~  W^y  +  Zkyy 


£Xy  =  Ufy  +  \^>X  +  W,XW,y  +  ZlCXy 


Where  kj  are  the  plate  curvatures  given  by 


(16) 


k**  =  -w,xx 


kyy  =  -W,yy 


k*y  —  ”2\y,xy 


(17) 


The  constitutive  relations  in  terms  of  the  generalized  strains  and 
stress  resultants  are  the  well  known  relations 


n  \  _  [A  Bl  fe°\ 
Ml  -  Lb  dJ  |k  j 


(18) 
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Considering  a  symmetric  laminate,  for  simplicity,  and  introducing 
Equation  (18)  into  Equation  (15)  and  B  =  0 


Du  w 

,xxxx  +  4Di6  w 

,xxxy  "*■  2(Dl2  +  2D66)  W>XXyy  +  4D26W>Xyyy 

+  4D26  W,Xyyy  +  D22  W.yyyy  =  NXW>XX  +  2NXy  W,Xy  + 

NyW.yy  +  P  (19) 

Next,  we  will  use  force  function  such  that  the  in-plane  equilibrium 
equations  are  satisfied 

NX=\|/,yy  Ny  =  \|TXX  Nxy  =  -\j/iXy  (20) 

Writing  the  compatibility  equation  in  terms  of  the  mid-surface 
strains 


Exx.yy  +  Eyy.xx  -  Exy.xy  =  W?xy  -  W,xx  W,yy  (21) 

Introducing  (20)  into  (18)  and  substituting  the  reference  surface 

pO  t 

strains  >i  into  (21) 


A22  V  ,xxxx  2A26  V  ,xxxy  +  (2A12  +  A55)  V,xxyy  ‘  2Ai6  \\f  ,xyyy  + 
A'n  V.yyyy  =  w?xy  -  w,xxw,yy 


(22) 


Equations  (19)  and  (22)  define  a  set  of  coupled  partial  differential 
equations.  These  are  the  postbuckling  nonlinear  equations  of  thin 
laminated  plates. 
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Fig. 3. -  Transverse  and  axial  stress-strain  curves,  ref.  [14]. 
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Fig.4  -  Convergence  study  of  postbuckling  response  of  clamped  plate  in 
the  form  of  plate  end-shortening. 


Fig. 5  -  Effect  of  laminate  thickness  on  plate  postbuckling  behavior  for 
linear  and  nonlinear  fiber  composite  materials.  End-shortening  versus 
axial  load  (clamped  plate,  l.%  imperfection,  Ny=0). 


Fig. 6  -  End-shortening  of  laminates  with  different  aspect-ratios  for 
linear  and  nonlinear  materials. (Clamped  plate,  \.%  imperfection, 
Ny/Nx=0.5). 
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Fig. 7  -  Effect  of  laminate  thickness  on  plate  postbuckhng  behavior  foi 
linear  and  nonlinear  fiber  composite  materials.  End-shortening  vet  sm 
axial  load.  (Clamped  plate,  l.%  imperfection,  Ny/Nx=l.). 
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Fig. 8  -  Effect  of  laminate  thickness  on  plate  postbuckling  behavior  loi 
linear  and  nonlinear  fiber  composite  materials.  End-shortening  versus 
axial  load.  (Clamped  plate,  l.%  imperfection,  Ny/Nx=2). 
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Fig. 9  -  End-shortening  of  thick-section  laminated  plate  for  linear  and 
nonlinear  material  models.  (l.%  imperfection,  Ny=0,  t/L=0.08). 
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Fig.  10  -  End -shortening  of  thick  section  laminated  plate  for  linear  and 
nonlinear  material  models.  (l.%  imperfection,  Ny/Nx-0.5,t/L=0.08). 
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Fig. 11  -  End-shortening  of  thick-section  laminated  plate  foi  linear  and 
nonlinear  material  models.  (l.%  imperfection,  Ny/Nx=2.,  t/L-O.UX) 
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Fig.  12  -  Effect  of  stress -biaxiality  ratio  on  postbuckling  response  for 
linear  and  nonlinear  material  models.  (l.%  imperfection,  t/L=0.08). 
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Fig.  14  -  Postbuckling  load  and  plate  central  deflection  for  thick-section 
and  different  material  models.  (l.%  imperfection,  Ny/Nx=0.5,  t/L=0.08). 


Fig.  16  -  End-shortening  of  thick-section  multi-layered  plate  for  linear 
and  nonlinear  material  model  (Clamped  plate,  l.%  imperfection,  Ny/Nx- 


Fig.  17  -  Lateral  deflection  of  multi-layered  thick  plate  for  different 
material  models.  (Clamped  plate,  l.%  imperfection,  Ny/Nx=l.,  t/L=0.08). 
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Fig.  18  -  Transverse  shear  behavior  in  a  clamped  [0/90/±45]izs  plate 
under  biaxial  compression. 
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Fig.  19  -  Transverse  shear  behavior  in  a  clamped  [0/90/±45]i2s  plate 
under  biaxial  compression. 


